Hydrodynamic modeling of granular flows in a modified Couette cell 
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We present simulations of granular flows in a modified Couette cell, using a continuum model 
recently proposed for dense granular flows. Based on a friction coefficient, which depends on an 
inertial number, the model captures the positions of the wide shear bands. We show that a smooth 
transition in velocity-profile shape occurs when increasing the height of the granular material, leading 
to a differential rotation of the central part close to the surface. The numerical predictions are in 
qualitative agreement with previous experimental results. The model provides predictions for the 
increase of the shear bands width when increasing the rotation rate. 

PACS numbers: 45.70.-n, 83. 50. Ax 



Divided materials such as granular material often ex- 
hibit a localization of deformations when slowly sheared 
[U, [2, 0| . This behavior produces shear bands that are 
generally narrow, typically 10 grain-diameters wide for a 
granular flow. They often take place close to boundaries 
or at the interface between a flowing and a static part. 
Quasi-static motions are experimentally observed inside 
these bands and can lead to different velocity proflles like 
exponential or Gaussian, depending on the system @, 
However a theoretical framework is still lacking to un- 
derstand them. A modified Couette configuration has 
recently brought the possibility of observing wide shear 
bands in granular media J5] and allows to more deeply 
test recent theories [H, Q • In this configuration sketched 
in Fig.[IJa), the bottom is split into a rotating disc (at a 
given angular velocity fio) and a ring fixed to the wall 
at rest. The three control parameters are the radius 
of the disc Rg, the filling thickness H and the angular 
velocity f2o- The experimental results are summarized 
here. For a small thickness compared to the disc ra- 
dius, the disc drives the rotation of a column up to the 
free surface as shown in Fig. mb). However, when the 
height is large, only a dome rotates deeply below the 
surface [Fig. [TJd)]. The transition occurs around the ra- 
tio H/Rs = 0.7 H, i, H, where two shear bands exist 
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[Fig. HJc)]. In this regime, a differential rotation (called 
precession in the following) appears: the central upper 
part rotates more slowly than the bottom [§, [l^]. A 
theoretical model, which is based on a principle of min- 
imization of the Coulomb friction dissipation, predicts 
infinitely sharp shear band and a hysteretic transition be- 
tween open and closed regimes [6|. Torok et al. recently 
introduce a random disorder, they find a smooth tran- 
sition and wide shear bands with corrects scaling laws 
when averaging over an ensemble of sharp bands 11]. 

In this paper, we present an alternative approach. A 
constitutive law has been proposed, which quantitatively 
describes some dense surface flows and instabilities [jjl, 
[Tsl . IT3 |. Could this hydrodynamic model describe the 
flow in the modified Couette cell? What are the shapes 
of bands and is there a smooth transition? 
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FIG. 1: a) Modified Couette cell: the bottom is split into a 
rotating disc and an outer ring at rest. From previous studies, 
in the stationary state, three regimes are identified depending 
on the filling height (the hatched areas rotate as quasi-solid 
bodies): (b) the open regime, (c) the intermediate regime 
when H/Rs ~ 0.7 and (d) the closed system. 



In the local model proposed in [l3] . the granular me- 
dia is described as an incompressible fluid governed by a 
visco-plastic constitutive law. Assuming an ensemble of 
grains of diameter d and density p^, the shear stress tXij 
is related to the strain rate jij as follows: 



with 



m(/)p 



I7l 



7y= 777y, (1) 



where P is the isotropic pressure, /i(/) is a friction 
coefficient, which depends on the inertial number / 
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the norm \X\ corresponds to y/XijXij /2 



and 7] is the effective viscosity. The friction coefficient 
starts from a minimum value fig and increases with / 
(see Ref. [l2| for details about the rheology). Although 
this model was not designed to reproduce quasi-static 
creeping flows (the constitutive law ensures the existence 
of a flnite yield stress |r| > figP), we expect it could 
contribute to understand what kind of phenomenon un- 
derlies the properties of shear bands in this geometry. We 
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assume the flow to be axisymmetric and compute the ax- 
ial velocity ug(r, z) = rQ(r, z), only in a radial section of 
the cell, governed by the following equations: 
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The simulations are made on a fixed grid (72 x 72). Us- 
ing a finite-difference scheme, we integrate the Cauchy 
equation for f2 from the material at rest until a steady 
state is reached with the following boundary conditions: 
a stress-free condition at the free surface and a no-slip- 
velocity condition at the wall and at the bottom (the an- 
gular velocity is then given by f2|(z=o) — for r < Rs 
and for r > Rs [i3|)- The set of parameters used in sim- 
ulations is similar to the experimental one found in Ref. 
[lol |. To simplify computations, we assume moreover the 
free surface to remain flat and horizontal: for low angular 
velocities the elevation is negligible and would reach 5% 
of the radius for a viscous fluid when f2o is multiplied by 
20. The pressure is thus the sum of a hydrostatic part 
and an inertial part due to the rotation. 

Figures [2ja-c) show typical angular- velocity profiles in 
the radial section of the cylinder for three different fill- 
ing heights. For H/d = 105 {H/Rs = 0.59), the central 
part, which rotates, emerges at the surface of the gran- 
ular material [Fig. [2i;a)] . For H/d = 130 {H/Rs = 0.73) 
in Fig.[2Kc), the rotating zone is reduced to a cupola near 
the rotating disc. Close to the transition [Fig. [^^b)], we 
observe two shear bands {H/Rs = 0.706): three quasi- 
rigid bodies are separated by narrow shear bands. The 
angular velocity of the upper central part is in-between 
and Oq. a first conclusion of this study is that the con- 
tinuum model qualitatively reproduces the experimen- 
tally observed phenomenon. Can it also describe the 
shape and the width of these shear bands? To study 
the shape, we define the position of shear bands in the 
bulk by the points where the angular velocity equals the 
average value of neighboring zones. Figure [Hd) shows 
the positions of the shear bands Zc{rc) for different filling 
heights. The continuous lines correspond to the extreme 
regimes. For increasing heights, the open shear band 
is more curved and is closer to the center. When the 
transition occurs, its radius at the free surface seems to 
reach a limit and the open shear band vanishes while the 
closed one appears below the surface. Both coexisting 
shear bands are shown by the dashed lines in Fig. [S^d). 
The dotted line corresponds to the empirical scaling of 
the shear position at the free surface (in the open state) 
from experiments of Fcnistein et al. @. It is in good 
agreement with the curve tips of our simulations, but 
a slight departure can be observed as obtained by the 
linger et aVs model 

Let us consider in more details the intermediate 
regime, where our simulations reveal the existence of two 
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FIG. 2; a-c) Stationary angular velocity profiles in a radial 
section of the Couette cell for Rs = 180d [Q.' = Q./lg/df^'^). 
a.) H = 105d {H/Rs = 0.59), the rotation reaches the free 
surface (open regime), h) H ^ 126. 5d {H/Rs = 0.706), 
this case corresponds to the transitory regime, c) H = ISOd 
{H/Rs — 0.73), above the transition we observe the closed 
regime, d) Position of the shear bands for different aspect ra- 
tio. Solid lines: H/Rs = 0.39, 0.5, 0.59, 0.63, 0.66 (open 
systems), 0.73, 0.78, 0.84 (closed systems). Dashed lines: 
H/Rs =0.700, 0.706, 0.714 (intermediate regime). Dotted 
line: experimental position of the shear band at the free sur- 
face from Ref. HI. 



shear zones bordering three rigid parts: the bottom ro- 
tates at rip, the lateral part is at rest and, in between, the 
upper central part rotates at a smaller angular velocity. 
To analyze the evolution of this phenomenon and to com- 
pare it with the experimental one , we plot in Fig. [3] 
the normalized difference between the angular velocity 
at the bottom and at the free surface on the symmetric 
axis: ujp = 1 — ^/^o\z=H,r=o- For a given angular veloc- 
ity, we obtain first that the transition between open and 
closed systems is continuous as observed in experiments 
[lot and DEM simulations 0, We do not find any 
hysteretic phenomenon unlike the predictions of Unger's 
model et al. Q. However, the transition, which occurs 
around H/R^ ~ 0.7, is sharper in our simulations than in 
experiments, although the simulated curve (squares) and 
the Fenistein et aVs experiments (crosses) correspond to 
the same set of parameters. This could be due to the 
inability of our model to capture the quasi-static proper- 
ties due to non-local effects (l6j . The important point is 
thus that we can reproduce the smooth transition using 
a shear-rate-dependent friction coefficient. As a conse- 
quence, the observed precession is the result of the su- 
perposition of two shear bands. When the rotation rate 
is increased by a factor 20 (fio = 3.2 rad.s"^), the tran- 
sition zone broadens as shown by the precession ratio ujp 
in Fig. [3] (triangles). We can notice moreover that, if the 
radial acceleration is negligible (left part of the equation 



3 



(On 



1 






A 


16 












0.8 








0.12 






X 






0.6 




X 




w 




A 

X 
















0.4 














A 




0.04 


0.2 
























A 









0.4 0.5 0.6 0.7 0.8 0.9 1 

H/Rs 

FIG. 3: Precession of the surface core: normalized difference 
of angular velocity cjp = 1 ~ Sl/f^o |z=H,r=o as a function of 
the rescaled height H/Rs- Comparison between simulations 
for two angular velocities (f2o = 0.16 rad.s"^ □ and Jlo ~ 
3.2 rad.s"^ A) and experimental results from Ref. ^To'] (x). 
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FIG. 4; Evolution of shear band widths as a function of the 
aspect ratio through the transition for two different angular 
velocities: a) Qo ~ 0.16 rad.s"^ o. b) fto = 3.2 rad.s"^ □. 
White (resp. black) symbols correspond to open (resp. 
closed) shear bands and the gray regions to the transition 
zone (0.03 < cup < 0.97). 



[3|) , equations governing the model are invariant when all 
lengths are normalized by and velocities are divided 

3/2 

by i?s . The latter exponent comes from the expres- 
sion of /. The rescaled variables are then H/Rs, W/Rg 

1 /2 

and il/Rs ■ The dependence on Qq can be translated in 
terms of Rg'. for very large systems, the extent of transi- 
tion should decrease for increasing R^ . This prediction is 
in agreement with Torok et aVs model 
diets an infinitely sharp transition for Rg 



which pre- 
oo. 

We now focus on the width of shear bands at the free 
surface and in the bulk, first for a low angular velocity 
(rio=0.16 rad.s^^). In the following, we will extract the 
band widths from the angular- velocity profiles in the sur- 
face and along the symmetric axis in the bulk fitted by 
an error function: fl/flo = 1/2-f l/2erf[(Xc-x)/W] (The 
axial shear band is very well fitted whereas the surface ve- 
locity profile is more linear). Following the dimensional 
analysis, we use the scaling W/Rg. The width of the 
open shear band, whose evolution with the filling height 
is shown in Fig. HJa), is very narrow. The width W in- 
creases with H / Rg (white circles) and suddenly decreases 
below the resolution of the grid close to the transition 
(gray region). After the transition, we plot the width 
of the vertical shear band along the axis of the cylinder. 
The vertical width in the closed system (black circles) 
is here larger and decreases with the height. Coinparing 
with experimental and numerical results [lol . [Tlj , we 
can notice that our open shear bands are very narrow 
but that the width of closed shear bands is of the same 
magnitude 

For a higher angular velocity (ilo=3.2 rad.s"^), in ad- 
dition to the enlargement of the transition range, both 
radial and vertical shear bands widen [Fig. IHb)]. The 
global behavior is similar to the one at low velocity, ex- 
cept that in that case, the variations of the widths within 
the transition state are more obvious: the widths of each 
band go to zero as they disappear. These variations of 

widths can be understood as follows: the rheology re- oc n%/R'^' ), implies instead that W/Rl^ is a function 



duces to a Coulomb friction law for small I (the normal- 
ized shear rate), so when the velocity difference between 
two neighbor zones decreases, the friction coefficient de- 
creases and the shear bands become narrower. These 
modifications of widths with respect to f2o are in con- 
trast with what is observed in direct simulations and in 
experiments where the velocity pattern is constant over 
a given range of rotation rate [a Q . Before concluding 
on this point, we study in more details the influence of 
the driven velocity on the width of the shear bands. 

In the following paragraph, we focus on two given as- 
pect ratios: an open shear band {H/Rg = 0.5) and a 
closed shear band {H/Rg — 0.84). First, we observe 
that the positions of shear bands are quite unchanged 
when increasing the rotation: the position of the open 
one moves very slightly toward the center. Second, we 
study the evolutions of widths with respect to the angu- 
lar velocity as shown in Fig. [51 For the open shear band, 
the width increases with a power law of angular velocity 
which exponent is a = 0.38. This relation means first 
that, in the limit of very low angular velocity, this shear 
band is infinitely narrow. This is consistent with the fact 
that our model reduces to a simple Coulomb friction law 
(a Drucker-Prager criterion) for I— s-0. In Unger et aVs 
model [fj], a constant friction law produces sharp shear 
bands. Concerning the closed shape, the widths have 
roughly the same behavior, but for low rotation rate we 
observe a slight departure from the power law. This be- 
havior could be attributed to numeric artifact that van- 
ishes when the size of the mesh grid decreases. Since 
this numerical effect does not alter the open shear band, 
we assume that the geometry could infiuence in different 
ways these two shear bands. From the dimensional anal- 
ysis, the scaling for the radial width reported in previous 

papers [W/R^^^ = /(H/Rg)] [13, d can not be derived: 
First we find that W/Rg is the relevant variable for the 
closed regime in the inertial regime, second the observed 
scaling law for the open shear band, W/Rg oc fig (i.e. 
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FIG. 5: Evolution of the width of a closed [H/Rs = 0.84 ■) 
and open {H/Rs = 0.50 o) shear band as a function of angular 
velocity Qq- The data are fitted by a power law of exponent 
0.38 (line). For the highest velocity and the open regime, the 
assumption of a flat free surface is no more valid. 

of H/Rs for a given ^q, with /3 = 1 - a/2 « 0.84. This 
discrepancy with previous results could come from the in- 
ability of the model to reproduce quasi-static phenomena, 
which govern the precise shape of shear bands. However, 
this dependance on the rotation rate might be seen ex- 
perimentally for higher angular velocities ((/) > 3 10"'^ 
for Qo > 0.02 ^/gjd). 

In conclusion, we have shown that a continuum model 
based on a hydrodynamic approach allows the computa- 
tion of the main properties of shear bands in the modi- 



fied Couette geometry: the threshold, which corresponds 
to a 3D generalization of the Coulomb friction thresh- 
old (|r| < figP), is sufficient to determine the position 
of shear bands, and the experimentally observed transi- 
tion in shapes is qualitatively reproduced. Second, we 
predict a smooth transition, which differs from the first- 
order transition predicted by Unger et aVs model. This 
comes from the increase of fi with the inertial number 
inside the shear band and from the superposition of two 
shear zones. Like in the plane-Couette geometry without 
gravity, the model predicts wide shear bands as soon as 
fi depends on /. However, the correct widths are not pre- 
dicted, specially the width of the radial shear band. This 
problem arises mainly from the quasi-static properties of 
the flow, but when we leave the quasi-static regime and 
enter in the dense regime (locally / w 0.05), some pre- 
dictions, which can be tested experimentally, are made 
concerning the evolution of these bands with respect to 
ilfl. In experiments at low shear, the effects of / may be 
hidden by other phenomenon. Another parameter, such 
as the one used in the Depken et al.'s shear-free-sheet 
model , may selects the properties of the velocity pro- 
files or the assumption of an isotropic pressure, which 
seems valid in rapid flows, may fail as it does in static 
piles. 

The author thanks O. Pouliquen, Y. Forterre for dis- 
cussions and carefully reading the manuscript and also 
M. van Hcckc for useful discussions. 
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